
# Table 1 : Panel A 

sampletab_1=dfbartik[dfbartik$testmonth=="May",]

ols_dist=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + avgpm25_9 + ", districtcovariates ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |0|leaid",sep="" ))
resultols_dist=felm(ols_dist,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)

ols_sort=as.formula(paste("cs_mn_all~ avgpm25_9 + cs_mn_all_lag_cohort + ",districtcovariates, "+ ", acs ," |factor(year) +factor(subject) + factor(grade)+factor(leaid) |0|leaid",sep="" ))
resultols_sort=felm(ols_sort,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)

# district + acs +weather  
ols_weathavg=as.formula(paste("cs_mn_all~ avgpm25_9 + cs_mn_all_lag_cohort + ", districtcovariates, " + ", acs ," + ", tempvarsavg ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |0|leaid",sep="" ))
resultols_weathavg=felm(ols_weathavg,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)

fileout=paste("output/Table1_PanelA.tex",sep="")

stargazer(resultols_dist,resultols_sort,resultols_weathavg, keep="avgpm25_9",digits=4, out=fileout)

############################# 
# Table 1 Panel B


#  district covariates
bartik_dist=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))

resultbartik_dist=felm(bartik_dist,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)


bartik_sort=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ",districtcovariates, "+ ", acs ," |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~",instname,")|leaid",sep="" ))
resultbartik_sort=felm(bartik_sort,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)


# district + acs +weather long 
bartik_weathavg=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates, " + ", acs ," + ", tempvarsavg ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))
resultbartik_weathavg=felm(bartik_weathavg,sampletab_1,weights=sampletab_1[[weightval]],na.action=na.omit)


# export 
fileout= paste("output/Table1_PanelB.tex",sep="")

stargazer(resultbartik_dist,resultbartik_sort,resultbartik_weathavg, keep="avgpm25_9",digits=4, add.lines=list(c("First-Stage F Stat",resultbartik_dist$stage1$iv1fstat$avgpm25_9[5], resultbartik_sort$stage1$iv1fstat$avgpm25_9[5], resultbartik_weathavg$stage1$iv1fstat$avgpm25_9[5]),
                                                                                                               c("Student Covariates", "Yes", "Yes", "Yes"),
                                                                                                               c("Local Economic Controls","No", "Yes", "Yes"),
                                                                                                               c("Weather Controls",   "No", "No", "Yes")
),
out=fileout )

##################################
# Table 1 Panel C : reading
sampletab_1c=sampletab_1[sampletab_1$subject=="rla",]

#  district covariates
bartik_dist=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates ,"  |factor(year)  + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))

resultbartik_dist=felm(bartik_dist,sampletab_1c,weights=sampletab_1c[[weightval]],na.action=na.omit)


bartik_sort=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ",districtcovariates, "+ ", acs ," |factor(year) + factor(grade)+factor(leaid) |(avgpm25_9 ~",instname,")|leaid",sep="" ))
resultbartik_sort=felm(bartik_sort,sampletab_1c,weights=sampletab_1c[[weightval]],na.action=na.omit)


# district + acs +weather long 
bartik_weathavg=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates, " + ", acs ," + ", tempvarsavg ,"  |factor(year)  + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))
resultbartik_weathavg=felm(bartik_weathavg,sampletab_1c,weights=sampletab_1c[[weightval]],na.action=na.omit)


# export 
fileout= paste("output/Table1_PanelC.tex",sep="")

stargazer(resultbartik_dist,resultbartik_sort,resultbartik_weathavg, keep="avgpm25_9",digits=4, add.lines=list(c("First-Stage F Stat",resultbartik_dist$stage1$iv1fstat$avgpm25_9[5], resultbartik_sort$stage1$iv1fstat$avgpm25_9[5], resultbartik_weathavg$stage1$iv1fstat$avgpm25_9[5]),
                                                                                                               c("Student Covariates", "Yes", "Yes", "Yes"),
                                                                                                               c("Local Economic Controls","No", "Yes", "Yes"),
                                                                                                               c("Weather Controls",   "No", "No", "Yes")
),
out=fileout )


##################################
# Table 1 Panel D : math
sampletab_1d=sampletab_1[sampletab_1$subject=="mth",]

#  district covariates
bartik_dist=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))

resultbartik_dist=felm(bartik_dist,sampletab_1d,weights=sampletab_1d[[weightval]],na.action=na.omit)


bartik_sort=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ",districtcovariates, "+ ", acs ," |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~",instname,")|leaid",sep="" ))
resultbartik_sort=felm(bartik_sort,sampletab_1d,weights=sampletab_1d[[weightval]],na.action=na.omit)


# district + acs +weather long 
bartik_weathavg=as.formula(paste("cs_mn_all~ cs_mn_all_lag_cohort + ", districtcovariates, " + ", acs ," + ", tempvarsavg ,"  |factor(year) +factor(subject) + factor(grade)+factor(leaid) |(avgpm25_9 ~ ",instname,")|leaid",sep="" ))
resultbartik_weathavg=felm(bartik_weathavg,sampletab_1d,weights=sampletab_1d[[weightval]],na.action=na.omit)


# export 
fileout= paste("output/Table1_PanelD.tex",sep="")

stargazer(resultbartik_dist,resultbartik_sort,resultbartik_weathavg, keep="avgpm25_9",digits=4, add.lines=list(c("First-Stage F Stat",resultbartik_dist$stage1$iv1fstat$avgpm25_9[5], resultbartik_sort$stage1$iv1fstat$avgpm25_9[5], resultbartik_weathavg$stage1$iv1fstat$avgpm25_9[5]),
                                                                                                               c("Student Covariates", "Yes", "Yes", "Yes"),
                                                                                                               c("Local Economic Controls","No", "Yes", "Yes"),
                                                                                                               c("Weather Controls",   "No", "No", "Yes")
),
out=fileout )

